Inferring Cellular Developmental Time

“From Multivariate to Longitudinal Data”

April 11, 2017

• Motivation
  • Single cell RNA-Seq
• Model Dataset
  • EDA
• Methods of estimating developmental time
  • PCA
  • Probablistic PCA; Bayesian PCA
  • Gaussian Process Latent Variable Modeling

• What happens in a cell such that it becomes a brain, toe, or a heart?
  • When do these decisions get made? “Who” makes them?
    
    
• How do the developmental trajectories of disease (leukemia / schizophrenia) differ from healthy individuals?
  
  
• Can we identify important transition points and the genetic signature underlying them?

\( \forall \) gene \( g \), fit OLS Regression with known timepoint \( t \) per cell–

\[ \log_2(g +1) \sim \beta_0 + \beta_1 t \]

\( \forall \) gene \( g \), fit regression with permuted timepoint \( t^* \) per cell

\[ \log_2(g +1) \sim \beta_0 + \beta_1 t^* \]

\( \forall \) gene \( g \), fit regression with permuted factor timepoint \( t^{**} \)

\[ \log_2(g +1) \sim \beta_0 + \beta_1 t^{**} \]

Given a matrix \( \mathbf{Y} \) of \( D \) genes (features) by \( n \) samples, determine a latent vector \( \mathbf{P} \) with dimension \( 1 \) x \( n \) that reflects the developmental trajectory of the \( n \) cells from the variance in \( D \) genes.

• \( D \) can be thought of has a higher dimension space, and we want to infer \( d \) (\( d \) < \( D \)) latent variables in the gene data.
  • One of the \( d \) latent variables ideally reflects developmental ordering.

\( \mathbf{Y} \) is a \( D \) x \( n \) matrix. Compute the covariance matrix–

\[ \mathbf{\Sigma} = E(\mathbf{YY^{T}}) - \mathbf{\mu \mu^T} \] where \( \mathbf{\mu} = E(\mathbf{Y}) \)

Then compute the spectral decomposition of \( \mathbf{\Sigma} \)

\[ \mathbf{\Sigma} a_j = \lambda_j a_j \] for \( j \in (1, ..., D) \)

Then \( a_j \) represent the eigenvectors of the data matrix \( \mathbf{Y} \).
Note the ordering of \( j \) is meaningful– \[ \lambda_1 \geq \lambda_2 \geq ... \geq \lambda_D \geq 0 \]

PCA by itself isn't satisfactory…



Ideas for improvements?

1) Quantifiying uncertainty

2) Non-linear latent variable inference

• Next several images taken from slides via Guy Wolf (Yale)

• These slides can be found here

\( \mathbf{Y} \) is a \( D \) x \( n \) matrix. \( \mathbf{x} \) is a \( d \) x \( n \) matrix (\( d < D \)). Want to relate \( \mathbf{x} \) and \( \mathbf{Y} \) and assume that the relationship is linear–

\[ \mathbf{Y = Wx} + \epsilon \]
where \( W \) is a \( D \) x \( d \) matrix.

By convention (after locating the matrix),

\[ \mathbf{x} \sim \mathcal{N} (0, \mathbf{I}) \]
where \( \mathbf{I} \) is the identity matrix of dimension \( d \) x \( d \).

Specify the error, \[ \epsilon \sim \mathcal{N} (0, \mathbf{\Psi}) \]
where we assume \( \mathbf{\Psi} \) to be a diagonal matrix \( D \) x \( D \).

• Assuming the diagonal structure of \( \mathbf{\Psi} \) means that all observertional correlations is due to the latent variables.
• In other words, \( Y_i \) for \( i \in (1, ..., n) \) are conditionally independent given \( \mathbf{x} \)

Then \( Y_i \) for \( i \in (1, ..., n) \),

\[ Y_i |x_i \sim \mathcal{N} (0, \mathbf{\Psi}) \]
Then we can integrate out the latent variables–

\[ p( Y_i |\mathbf{W}, \mathbf{\Psi})= \int p( Y_i | x_i , \mathbf{W, \Psi}) p(x_i) dx_i \]

Under this specification, we can write the MVN distribution for our observations–

\[ \mathbf{Y} \sim \mathcal{N} (0, \mathbf{C}), \mathbf{C} = \mathbf{WW^T }+ \mathbf{\Psi} \]

Young-Whittle Factor Analysis (1950s)

\( \psi_i \) element of the diagonal of \( \Psi \); add constraint \( \psi_i = \sigma^2 \)

Assume \( \sigma^2 \) known, MLE yields same \( W \) as PCA

\[ \mathcal{L} = \tfrac{-N}{2}\{ D\log(2\pi) + \log|\mathbf{C}| + \textrm{tr}(\mathbf{C}^{-1}\mathbf{\Sigma}) \} \]
\[ \mathbf{C} = \mathbf{WW^T }+ \mathbf{\Psi} \]
\[ \mathbf{x} = \mathbf{C^{-1}W^{T}Y} \]

• Assuming \( \sigma^2 \) may not be reasonable; want to estimate it from the data (keep in likelihood)
• Can estimate \( \mathbf{W}, \sigma^2 \) using EM & with a prior over \( \mathbf{x} \)

For \( \mathbf{W_{\textrm{MLE}}} \),

\[ \sigma_{MLE}^2 = \dfrac{1}{D-d} \sum_{j = d +1}^D \lambda_j \]

• Similarly, we can integrate over \( \mathbf{W} \) given a prior, yielding
  
  \[ \mathbf{x} \sim \mathcal{N} (\mathbf{C^{-1}W^{T}Y}, \sigma^2C^{-1}) \]
  

• \( \alpha_i \) represents the inverse of the variance
  • for large \( \alpha_i \), contribution of latent factor \( i \) is low.
• Solved by EM / similar algorithm

\[ \mathbf{Y} \sim \mathcal{N} (0, \mathbf{C}), \mathbf{C} = \mathbf{WW^T }+ \mathbf{\Psi} \]
\[ \mathbf{C}(\mathbf{W_i},\mathbf{W_j}) = \mathbf{W_i^T}\mathbf{W_j} + \sigma^2\delta_{ij} \]
\[ \mathbf{C}(\mathbf{W_i},\mathbf{W_j}) = \theta_{rbf}\exp\Big(\dfrac{-\gamma}{2}(\mathbf{W_i}-\mathbf{W_j})^{T}(\mathbf{W_i}-\mathbf{W_j}) \Big) + ... \]

• (2) being a special case (linear, iid) of (3)
• Computationally more challenging, but there are fast algorithms out there

• GPLVM provide both a probablistic and non-linear latent variable inference structure
  
  • All other published algorithms provide point estimates, difficult to ascertain uncertainty
    
    
• Under this framework, Bayesian priors are allowed, which have been useful in published studies
  
  
• Tools and methods useful in many data analyses contexts (including machine learning)
  • Motivated here by single cell analysis